capture log close
clear
capture program drop _all
set more off
timer clear 1
timer on 1
set seed 123

**************** 
***  PATHS   ***  
****************
*BCCR: 1 if at Central Bank
global BCCR=1
if $BCCR == 1 {
	global path "Y:/Projects/MNC_Suppliers/Suppliers_Alfaro-Urena_Manelici_Vasquez"
}

if $BCCR == 0 {
	global path "~/Documents/GitHub/Suppliers_Alfaro-Urena_Manelici_Vasquez"
}
cd $path 

*log file
local logfile "${path0}/log_files/6-appendix_c"
set linesize 255
log using `logfile', replace
set linesize 255

/*
RESULTS IN APPENDIX C OF THE PAPER
*/

********************** 
***  INPUT FILES   ***  
**********************
global analysis_data ""processed_data/analysis_data.dta""
global mnc_countries ""raw_data/mnc_countries_sample.dta""
global management ""raw_data/wmsdata_2004_2015.xlsx""
global country_codes ""raw_data/country_codes.xlsx""
global revec_group ""raw_data/revec_group.dta""
global transactions ""raw_data/trans_clean_corresp.dta""
global corresp ""raw_data/2018-08-14 Correspondencias de n° de identificación.xlsx""
global revec_corresp ""raw_data/revec_corresp.dta""
global imports_customs ""raw_data/imports_data.dta""
global exports_customs ""raw_data/exports_data.dta""
global MNCs ""processed_data/MNC_sample.dta"" 
global mnc_countries ""raw_data/mnc_countries_sample.dta""
global events_data ""processed_data/events_data.dta""
*dummies to show in tables and graphs 	
global dummies "D_m4 D_m3 D_m2 D_0 D_1 D_2 D_3 D_4"

********************************************************************************
* main program: outline of the full code
********************************************************************************
prog main

* Table C1
* Domestic Firms (Weakly) Reduce their Markups after Starting to Supply to MNCs
		markups

* Figure C1
* Split of Corporate Sales by Old or New, Domestic or MNC Buyers
		sh_to_other_MNCs

* Figure C2 
* Decomposition of Total Sales for First-time Suppliers to MNCs
		area_graph_sales

* Tables C2 and C3
* Heterogeneity by Number of MNC Events Experienced
		multiple_events

		
* Table C4
* Heterogeneity in TFP Gains Based on the Characteristics of the First MNC Buyer
		tfp_mnc_char_hete

* Table C5
* Heterogeneity in TFP Gains Based on the HQ Country of the First MNC Buyer
		het_mnc_countries

end


********************************************************************************
*Table C1: Domestic Firms (Weakly) Reduce their Markups after Starting to Supply to MNCs
********************************************************************************
{
prog markups 
display "Table C1: markups estimation and event time"

*DLW GMM BASED ON DLW (2012) CODE
	dlw_gmm

*MATA GMM routine DWL
	quietly do ${path0}/all_codes/Z-GMM_mata.do
	
*markups table DLW	
	markups_event

*markups Sampi et al. Data prep before matlab code
	main_sampi
	
end

******
*DLW (2012) GMM
******
prog dlw_gmm
display "markups: dlw gmm" 
quiet{

use $analysis_data , clear

clear mata
gen temp_l=l
gen temp_k=k
gen temp_m=m
drop l* k* m*
rename temp_l l
rename temp_k k
rename temp_m m 

* higher order terms on inputs
local M=3
local N=3
forvalues i=1/`M' {
gen l`i'=l^(`i')
gen m`i'=m^(`i')
gen k`i'=k^(`i')
*interaction terms
forvalues j=1/`N' {
gen l`i'm`j'=l^(`i')*m^(`j')
gen l`i'k`j'=l^(`i')*k^(`j')
gen k`i'm`j'=k^(`i')*m^(`j')
}
}
gen lkm=l*k*m
gen l2k2m2=l2*k2*m2
gen l3k3m3=l3*k3*m3
* generate interaction terms of all (l,m,k) 
*with exportdummy:call all terms e*(.) below
*------------------------------------------------------------------------------*
* INITIAL VALUES
* CD
quietly reghdfe va m l k, absorb(IDs year)
scalar blols=_b[l]
scalar bkols=_b[k]
scalar bmols=_b[m]

*------FIRST STAGE  -----------------------------------------------------------*
*xi: reg y l* m* k* i.year
quietly reghdfe va m* l* k*, absorb(year IDs) res(epsilon)
predict phi
label var phi "phi_it 
label var epsilon "measurement error first stage
xtset IDs year, yearly
gen phi_lag=L.phi
*------------------------------------------------------------------------------*
gen l_lag=L.l
gen k_lag=L.k
gen m_lag=L.m
gen l_lag2=l_lag^2
gen k_lag2=k_lag^2
gen m_lag2=m_lag^2
gen l_lagk_lag=l_lag*k_lag
gen l_lagm_lag=l_lag*m_lag
gen k_lagm_lag=k_lag*m_lag
gen lm=l*m
gen lk=l*k
gen km=k*m
gen l_lagk=l_lag*k
gen l_lagm=l_lag*m
*---COMPUTE CORRECTED SHARES---------------------------------------------------*
gen y_c=va-epsilon
gen va_c=exp(y_c)
gen alpha_l=salarios/va_c/1000
*------------------------------------------------------------------------------*
sort IDs year
gen const=1
drop if y==.
drop if l_lag==.
drop if k==.
drop if m==.
drop if phi==.
drop if phi_lag==.

}
end

******
*saving matrix of coefficients
******
program dlw, rclass
preserve 
sort IDs year
mata DLW()
return scalar beta_l=abs(beta_dlw[1,2])
return scalar beta_k=abs(beta_dlw[1,3])

gen mu_= beta_dlw[1,2]/alpha_l
quietly summ mu_, det
local per1=`r(p1)'
local per99=`r(p99)'
keep if mu_<`per99' & mu_>`per1'
quietly summ mu_, det
return scalar mu=`r(mean)'
return scalar obs_=_N
end

******
*saving matrix of coefficients
******
prog markups_event
display "markups: table" 
quiet{
*program that saves coefficients
dlw
gen mu_= beta_dlw[1,2]/alpha_l
quietly summ mu, det
local per1=`r(p1)'
local per99=`r(p99)'
replace mu=. if mu>=`per99' | mu<=`per1'

ds D_m4-D_4
global dummies `r(varlist)'

***  TABLE: MARKUPS GRAPH  
estimates clear
capture drop red
* COLUMN 1: NEVER MATCHED + EVENTUALLY MATCHED, VAR: sales
global REG "reghdfe mu  D_* if past!=1 , absorb(year#sector4#province IDs) vce(cluster province#sector_g)"
$REG resid(red)
distinct IDs if red!=.
scalar n_firms=r(ndistinct)
capture drop dep
* exchange rate around 500 colones per US dollar
gen dep=mu
quietly summ dep if red!=.
scalar mean_dep=r(mean)
scalar sd_dep=r(sd)
estadd scalar mean_dep
estadd scalar sd_dep
scalar n_obs=e(N)
scalar n_fe=e(df_a_initial)
display `=n_fe' 
estadd scalar n_firms 
estadd scalar n_obs
estadd scalar n_fe
estadd local Never_Matched "YES"
addfe Never_Matched
estimates store m_1
capture drop red
* COLUMN 2: EVENTUALLY MATCHED, VAR: sales
global REG "reghdfe mu  D_*  if event!=. , absorb(year#sector4#province IDs) vce(cluster province#event)"
$REG resid(red)
distinct IDs if red!=.
scalar n_firms=r(ndistinct)
capture drop dep
gen dep=mu
quietly summ dep if red!=.
scalar mean_dep=r(mean)
scalar sd_dep=r(sd)
estadd scalar mean_dep
estadd scalar sd_dep
scalar n_obs=e(N)
scalar n_fe=e(df_a_initial)
display `=n_fe' 
estadd scalar n_firms 
estadd scalar n_obs
estadd scalar n_fe
estadd local Never_Matched "NO"
estimates store m_2
capture drop red
* Prepare estimates for -estout-
	estfe m_*, labels(IDs "Firm FE" year#sector4#province "Year-4DSect-Prov FE")
	return list	
}	
esttab m_* using ${path0}/results/3-appendix_c/Supplementary_Table_C1a.tex, ///
 indicate(`r(indicate_fe)' Never_Matched)  star(* 0.10 ** 0.05 *** 0.01) ///
 varwidth(23) keep($dummies ) ///
 cells("b(fmt(3) star )" "se(par fmt(3))" ) ///
 stats(r2_a mean_dep sd_dep n_obs n_fe n_firms , fmt(a2) labels("Adjusted R$^2$" "Mean Dep. Var. (level)" "SD Dep. Var. (level)" "\# Observations" "\# Fixed Effects" "\# Firms" ))   ///
 nodepvars nomtitles alignment(c) replace noconstant label nodepvar collabels(none) compress type
end

********************************************************************************
* main program: outline of the full code
********************************************************************************
prog main_sampi

*Impact on markups by Sampi: preparing inputs for matlab
*Prepares firms' yearly data (labor share, total sales, nominal stock of capital and number of employees) that enters as input in Sampi's et al. (2021) MATLAB code. It produces the inputs: `matlab_input_${var}_${method}.csv`. 

*Sampi method (to estimate Table_C1_b_pooled)
	markups_sampi_global
*Sampi method, dropping never treated and using eventually treated as controls (to estimate Table_C1_c_pooled)
	markups_sampi_global_never		
*Sampi method (relative time) (to estimate Table_C1_b_yearly)
	markups_sampi_rel_med_aleat	
*Sampi method (relative time), dropping never treated and using eventually treated as controls (to estimate Table_C1_c_yearly)
	mps_smpi_rel_med_nver_aleat		
	
end


********************************************************************************
*Sampi method
********************************************************************************
prog markups_sampi_global

use $analysis_data , clear

*monetary from millions to colones; salarios from thousands to colones
global list_var "va_ ingresosir_ total_activo_neto_"
foreach var in $list_var {
global var "`var'"
replace $var = $var *1000000 
}
replace salarios_ = salarios_*1000

*filling missing time for each firm 
keep ID year event trabaj_ va_ ingresosir_ total_activo_neto_ salarios_
rename (trabaj_ va_ ingresosir_ total_activo_neto_ salarios_) (trabaj va ingresosir capital salarios)
destring ID, replace
tsset ID year
*filling missing time info
tsfill, full
tostring ID, replace


sort ID year
bys ID: egen temp_event = max(event)
drop event
rename temp_event event
*labor share
gen labor_sh_2 = salarios/va
replace labor_sh_2 = 1 if labor_sh_2>1 & labor_sh_2!=.
*event-regulation dummy: 1 if event active or after event, 0 otherwise
gen D =(year>=event & event!=.)
rename ingresosir ventas


*dealing with rows of NaNs: Fxt
global list_var "capital ventas labor_sh_2 D trabaj"
foreach var in $list_var {
global var "`var'"
drop if $var ==. 
}
global list_var "capital ventas labor_sh_2 trabaj"
foreach var in $list_var {
global var "`var'"
drop if $var ==0 
}

*math condition in sampi-matlab
drop if capital<trabaj & capital!=. & trabaj!=.

bys ID: gen count_yr = _N
drop if count_yr!=10


*long to wide
global list_var "capital ventas labor_sh_2 D trabaj"
foreach var in $list_var {
preserve
global var "`var'" 
keep ID year $var
rename $var var
reshape wide var, j(year) i(ID)
drop ID
export delimited using "temp/matlab_input_${var}_global.csv", replace
restore
}


end


********************************************************************************
*Sampi method, dropping never treated and using eventually treated as controls 
********************************************************************************
prog markups_sampi_global_never

use $analysis_data , clear

*monetary from millions to colones; salarios from thousands to colones
global list_var "va_ ingresosir_ total_activo_neto_"
foreach var in $list_var {
global var "`var'"
replace $var = $var *1000000 
}
replace salarios_ = salarios_*1000

*filling missing time for each firm 
keep ID year event trabaj_ va_ ingresosir_ total_activo_neto_ salarios_
rename (trabaj_ va_ ingresosir_ total_activo_neto_ salarios_) (trabaj va ingresosir capital salarios)
*dropping never treated and using eventually treated as controls
drop if event==.
gen event_control_treated=event
tempfile eventually
save `eventually', replace
*creating controls
replace ID = ID+ID
replace event=.
*appending treated
append using `eventually'
destring ID, replace
sort ID year

tsset ID year
*filling missing time info
tsfill, full
tostring ID, replace

sort ID year
bys ID: egen temp_event = max(event)
drop event
rename temp_event event
bys ID: egen temp_event_control_treated = max(event_control_treated)
drop event_control_treated
rename temp_event_control_treated event_control_treated
*labor share
gen labor_sh_2 = salarios/va
replace labor_sh_2 = 1 if labor_sh_2>1 & labor_sh_2!=.
*event-regulation dummy: 1 if event active or after event, 0 otherwise
gen D =(year>=event & event!=.)
rename ingresosir ventas


*dealing with rows of NaNs: Fxt
global list_var "capital ventas labor_sh_2 D trabaj"
foreach var in $list_var {
global var "`var'"
drop if $var ==. 
}
global list_var "capital ventas labor_sh_2 trabaj"
foreach var in $list_var {
global var "`var'"
drop if $var ==0 
}

*math condition in sampi-matlab
drop if capital<trabaj & capital!=. & trabaj!=.

bys ID: gen count_yr = _N
drop if count_yr!=10


*long to wide
global list_var "capital ventas labor_sh_2 D trabaj"
foreach var in $list_var {
preserve
global var "`var'" 
if "$var"!="D"{
replace $var =1.123456 if year>=event_control_treated & event==.	
}
keep ID year $var
rename $var var
reshape wide var, j(year) i(ID)
drop ID
export delimited using "temp/matlab_input_${var}_global_never.csv", replace
restore
}


end


********************************************************************************
*Sampi method (relative time)
********************************************************************************
prog markups_sampi_rel_med_aleat

set seed 123
use $analysis_data , clear

*monetary from millions to colones; salarios from thousands to colones
global list_var "va_ ingresosir_ total_activo_neto_"
foreach var in $list_var {
global var "`var'"
replace $var = $var *1000000 
}
replace salarios_ = salarios_*1000

*filling missing time for each firm 
keep ID year event trabaj_ va_ ingresosir_ total_activo_neto_ salarios_
rename (trabaj_ va_ ingresosir_ total_activo_neto_ salarios_) (trabaj va ingresosir capital salarios)
destring ID, replace
sort ID year

bys ID: egen min_yr = min(year)
bys ID: egen max_yr = max(year)
gen med_yr_temp = runiform(min_yr,max_yr)
replace med_yr_temp=round(med_yr_temp, 1)
replace med_yr_temp=. if year!=min_yr
bys ID: egen med_yr = max(med_yr_temp)
gen no_event=(event==.)
replace event = med_yr if event==.
drop min_yr max_yr med_yr_temp med_yr

tsset ID year
*filling missing time info
tsfill, full
tostring ID, replace


sort ID year
bys ID: egen temp_no_event = max(no_event)
drop no_event
rename temp_no_event no_event
bys ID: egen temp_event = max(event)
drop event
rename temp_event event
*relative time
gen time = year - event
*labor share
gen labor_sh_2 = salarios/va
replace labor_sh_2 = 1 if labor_sh_2>1 & labor_sh_2!=.
*event-regulation dummy: 1 if event active or after event, 0 otherwise
gen D =(time>=0 & no_event==0)
rename ingresosir ventas

*dealing with rows of NaNs: Fxt
global list_var "capital ventas labor_sh_2 D trabaj"
foreach var in $list_var {
global var "`var'"
drop if $var ==. 
}
global list_var "capital ventas labor_sh_2 trabaj"
foreach var in $list_var {
global var "`var'"
drop if $var ==0 
}

*math condition in sampi-matlab
drop if capital<trabaj & capital!=. & trabaj!=.

bys ID: gen count_yr = _N
drop if count_yr!=10


keep if time>=-4 & time<=4
tostring time, replace
replace time = "m4" if time=="-4"
replace time = "m3" if time=="-3"
replace time = "m2" if time=="-2"
replace time = "m1" if time=="-1"
replace time = "p0" if time=="0"
replace time = "p1" if time=="1"
replace time = "p2" if time=="2"
replace time = "p3" if time=="3"
replace time = "p4" if time=="4"


*long to wide
global list_var "capital ventas labor_sh_2 D trabaj"
foreach var in $list_var {
preserve
global var "`var'" 
keep ID time $var
rename $var var
reshape wide var, j(time) i(ID) string

if "$var"=="D"{	
global varlist "varm4 varm3 varm2 varm1 varp0 varp1 varp2 varp3 varp4"	
foreach x in $varlist{
  replace `x' = 0 if `x'==. 
}	
global varlist "varp1 varp2 varp3 varp4"	
foreach x in $varlist{
  replace `x' = 1 if varp0==1 
}
}

drop ID
order varm4 varm3 varm2 varm1 varp0 varp1 varp2 varp3 varp4 
export delimited using "temp/matlab_input_${var}_rel_med_aleat.csv", replace
restore
}



end


********************************************************************************
*Sampi method (relative time), dropping never treated and using eventually treated as controls 
********************************************************************************
prog mps_smpi_rel_med_nver_aleat

set seed 123
use $analysis_data , clear

*monetary from millions to colones; salarios from thousands to colones
global list_var "va_ ingresosir_ total_activo_neto_"
foreach var in $list_var {
global var "`var'"
replace $var = $var *1000000 
}
replace salarios_ = salarios_*1000

*filling missing time for each firm 
keep ID year event trabaj_ va_ ingresosir_ total_activo_neto_ salarios_
rename (trabaj_ va_ ingresosir_ total_activo_neto_ salarios_) (trabaj va ingresosir capital salarios)
*dropping never treated and using eventually treated as controls
drop if event==.
gen no_event = 0
tempfile eventually
save `eventually', replace
*creating controls
replace ID = ID+ID
replace no_event=1
*appending treated
append using `eventually'
destring ID, replace
sort ID year

bys ID: egen min_yr = min(year)
gen med_yr_temp = runiform(min_yr,event)
replace med_yr_temp=round(med_yr_temp, 1)
replace med_yr_temp=. if year!=min_yr
bys ID: egen med_yr = max(med_yr_temp)
gen event_control_treated = event
replace event=. if no_event==1
replace event = med_yr if event==.
drop min_yr med_yr_temp med_yr

tsset ID year
*filling missing time info
tsfill, full
tostring ID, replace


sort ID year
bys ID: egen temp_no_event = max(no_event)
drop no_event
rename temp_no_event no_event
bys ID: egen temp_event = max(event)
drop event
rename temp_event event
bys ID: egen temp_event_control_treated = max(event_control_treated)
drop event_control_treated
rename temp_event_control_treated event_control_treated
*relative time
gen time = year - event
*labor share
gen labor_sh_2 = salarios/va
replace labor_sh_2 = 1 if labor_sh_2>1 & labor_sh_2!=.
*event-regulation dummy: 1 if event active or after event, 0 otherwise
gen D =(time>=0 & no_event==0)
rename ingresosir ventas

*dealing with rows of NaNs: Fxt
global list_var "capital ventas labor_sh_2 D trabaj"
foreach var in $list_var {
global var "`var'"
drop if $var ==. 
}
global list_var "capital ventas labor_sh_2 trabaj"
foreach var in $list_var {
global var "`var'"
drop if $var ==0 
}

*math condition in sampi-matlab
drop if capital<trabaj & capital!=. & trabaj!=.

bys ID: gen count_yr = _N
drop if count_yr!=10

keep if time>=-4 & time<=4
tostring time, replace
replace time = "m4" if time=="-4"
replace time = "m3" if time=="-3"
replace time = "m2" if time=="-2"
replace time = "m1" if time=="-1"
replace time = "p0" if time=="0"
replace time = "p1" if time=="1"
replace time = "p2" if time=="2"
replace time = "p3" if time=="3"
replace time = "p4" if time=="4"


*long to wide
global list_var "capital ventas labor_sh_2 D trabaj"
foreach var in $list_var {
preserve
global var "`var'" 
replace $var =1.123456 if year>=event_control_treated & no_event==1
keep ID time $var
rename $var var
reshape wide var, j(time) i(ID) string

if "$var"=="D"{	
global varlist "varm4 varm3 varm2 varm1 varp0 varp1 varp2 varp3 varp4"	
foreach x in $varlist{
  replace `x' = 0 if `x'==.
  replace `x' = 0 if `x'==float(1.123456)
}	
global varlist "varp1 varp2 varp3 varp4"	
foreach x in $varlist{
  replace `x' = 1 if varp0==1 
}
}

drop ID
order varm4 varm3 varm2 varm1 varp0 varp1 varp2 varp3 varp4 
export delimited using "temp/matlab_input_${var}_rel_med_never_aleat.csv", replace
restore
}
end
}
********************************************************************************
* Figure C1
* Split of Corporate Sales by Old or New, Domestic or MNC Buyers
********************************************************************************
prog sh_to_other_MNCs
display "share to other buyers (old and new, domestic and MNCs)"
quiet{
*NEW/OLD SUPPLIERS, SECTORS, ETC
estimates clear
use $analysis_data, clear

*  	DUMMIES
*dummies to show in tables  	
ds D_m4-D_4
global dummies `r(varlist)'
*in the main sample 
reghdfe y D* if past!=1, absorb(year#province#sector4 IDs) resid(temp_res)
bys ID: egen res=mean(temp_res)
replace res=(res!=.)

keep event ingresos costo_de_ imports past province sector* ID year y res D* IDs

*sellers
preserve
keep if res==1
rename * s_*
rename s_ID seller
rename s_year year
tempfile as_seller
save `as_seller', replace
restore 

*buyers
drop D*
rename * b_*
rename b_ID buyer
rename b_year year
tempfile as_buyer
save `as_buyer', replace

* MNCS LIST
use $MNCs, clear
keep ID
rename ID buyer
tempfile MNC_buyer
save `MNC_buyer', replace
rename buyer seller
tempfile MNC_seller
save `MNC_seller', replace

* TRANSACTIONS
use $transactions, clear
drop if pot==1
drop pot
*REAL VALUE
local trans_vari "trans"
foreach v of local trans_vari{
forvalue i=2008/2017{
replace `v'=(100*(`v')/${ipc`i'})/1000000 if year==`i'
}
}

* MERGING BUYERS AND SELLERS
merge m:1 buyer year using `as_buyer'
bys buyer: egen buyer_data=max(_m)
replace buyer_data=(buyer_data==3)
drop _m
merge m:1 seller year using `as_seller'
bys seller: egen seller_data=max(_m)
replace seller_data=(seller_data==3)
drop _m
keep if buyer_data==1 | seller_data==1
* WHICH ARE MNCS
merge m:1 buyer using `MNC_buyer'
gen buyer_mnc=(_m==3)
drop if _m==2
drop _m

* buyers that triggered the event
gen trig=(year==s_event & buyer_mnc==1)
bys seller buyer: egen mnc_trig=max(trig)
gen tr_not_event=trans_*(1-mnc_trig)
gen l_trans_not_event=log(1000*tr_not_event)
egen sellers=group(seller)
egen buyers=group(buyer)
gen time= year-s_event

*old buyers
gen temp_old_b=(time<=0) 
bys buyer seller: egen old_b=max(temp_old_b)
drop temp*
scalar last_year=4
scalar first_year=-4
gen time_reg=time
replace time_reg=`=last_year' if time>=`=last_year' & time!=.
replace time_reg=`=first_year' if time<=`=first_year'
 }
***
*** figure with percentages depending on buyer types and table
***
quiet{
preserve
gen buyer_type_1= (old_b==1 & mnc_trig!=1)
gen buyer_type_2= (mnc_trig==1)
gen buyer_type_3= (old_b==0 & buyer_mnc!=1)
gen buyer_type_4= (buyer_mnc==1 & mnc_trig!=1)
assert buyer_type_1 + buyer_type_2 + buyer_type_3 + buyer_type_4 ==1
forvalue i=1/4{
gen n_buyer_type_`i'=buyer_type_`i'
replace buyer_type_`i'=buyer_type_`i'*trans_
}
*total trans
bys seller year: egen tot_trans=total(trans_)
collapse (sum) n_buyer_type_* buyer_type_* (mean) tot_trans s_ingresos , by(seller year time)
forvalue i=1/4{
gen sh_buyer_type_`i'=buyer_type_`i'/s_ingresos
replace sh_buyer_type_`i'=1 if sh_buyer_type_`i'>1 & sh_buyer_type_`i'!=.
}
gen time_reg=time
replace time_reg=. if time<0 | time>4
label def sa 0"0" 1"+1" 2"+2" 3"+3" 4"+4" 
label val time_reg sa
estpost tabstat  sh_buyer_type_2, statistics(count mean sd) by(time_reg) nototal  
est store share1
estpost tabstat  sh_buyer_type_4, statistics (count mean sd) by(time_reg) nototal  
est store share2
estpost tabstat  n_buyer_type_4, statistics(count mean sd) by(time_reg) nototal  
est store n_others
}

quiet{
collapse (sum) buyer_type_* tot_trans , by(time)
forvalue i=1/4{
replace buyer_type_`i'=buyer_type_`i'/tot_trans
}
keep if time>=-4 & time<=4

graph bar buyer_type*  , over(time ) stack percent ///
bar(1, color(navy) fintensity(inten100)) ///
bar(2, color(maroon) fintensity(inten100)) ///
bar(3, color(navy) fintensity(inten40)) ///
bar(4, color(maroon) fintensity(inten40)) ///
ytitle("Share of Total Sales") /// 
yscale(titlegap(*+10)) ///
b1title("Year Since Event") ///
graphregion(color(white)) bgcolor(white) /// 
legend(order(1 "Old DOM" 2 "First MNC" ///
3 "New DOM" 4 "Other MNCs") rows(1) region(lwidth(none)))    
graph export "${path0}/results/3-appendix_c/Supplementary_Figure_C1a.eps", as(eps) replace
restore

***
*** random event for those not treated
***
preserve 
bys seller: egen min_t=min(time)
drop if min_t!=.
drop if s_past==1
collapse year , by(seller)
gen event=ceil(2009+6*runiform())
keep seller event
tempfile other_sellers
save `other_sellers', replace
restore 

preserve
merge m:1 seller using `other_sellers'
keep if _m==3
capture drop time
gen time=year-event
*old buyers
capture drop old*
gen temp_old_b=(time<=0) 
bys buyer seller: egen old_b=max(temp_old_b)
*total trans
bys seller year: egen tot_trans=total(trans_)
*old and new buyers
gen buyer_type_1= (old_b==1)*trans_
gen buyer_type_2= (old_b!=1)*trans_
collapse (sum) buyer_type_* (mean) tot_trans , by(seller year time)
collapse (sum) buyer_type_* tot_trans , by(time)
forvalue i=1/2{
replace buyer_type_`i'=buyer_type_`i'/tot_trans
}
keep if time>=-4 & time<=4
graph bar buyer_type*  , over(time ) stack percent ///
bar(1, color(navy) fintensity(inten100)) ///
bar(2, color(navy) fintensity(inten40)) ///
ytitle("Share of Total Sales") /// 
b1title("Year Since Event") ///
graphregion(color(white)) bgcolor(white) /// 
legend(order(1 "Old DOM" ///
2 "New DOM") rows(1) region(lwidth(none)))  
graph export "${path0}/results/3-appendix_c/Supplementary_Figure_C1b.eps", as(eps) replace 
restore
}
end


********************************************************************************
* Figure C2
* Decomposition of Total Sales for First-time Suppliers to MNCs
********************************************************************************
prog area_graph_sales
display "figure C3: Sales decomposition"
quiet{
use $events_data , clear
rename ID seller
tempfile suppliers
keep seller ingresos event year exports
save `suppliers', replace

use $MNCs, clear
rename ID buyer
keep buyer
tempfile mnc
save `mnc', replace
* balance sheet 
use $revec_group , clear
replace ipubgrupo=1 if inolucrogrupo==1
keep ID foreign_list ipubgrupo firm_AE firm_ciiu nombre ZF
rename ID buyer
merge 1:1 buyer using `mnc'
gen MNC=(_m==3)
drop _m
replace foreign=1 if MNC==1
tempfile buyers
save `buyers', replace

* transactions 
use $transactions
drop if pot==1
merge m:1 seller year using `suppliers'
keep if _m==3
drop _m

bys seller: egen temp_event=mean(event)
drop event
rename temp_eve event

merge m:1 buyer using `buyers'
keep if _m==3
gen time=year-event
keep if time!=.
 
sort time trans_
by time : gen n_tot=_N
bys time : gen num=_n
replace num=num/n_tot
drop if num>0.995 & MNC!=1 
drop if num>0.99 & MNC!=1 & time==-2
drop num n_tot
***
*TRANSACTIONS TYPE
***
** total transactions
bys seller year: egen tot_trans=total(trans_)

** trans with domestic
gen temp =trans_*(1-forei)*(1-ipub)
bys seller year: egen tot_dom=total(temp)
drop temp

** trans with gov
gen temp =trans_*ipub
bys seller year: egen tot_gov=total(temp)
drop temp

** trans with MNC
gen temp =trans_*MNC
bys seller year: egen tot_MNC=total(temp)
drop temp

** trans with fore non mnc
gen temp =trans_*foreign_*(1-MNC)
bys seller year: egen tot_fore_nonmnc=total(temp)
drop temp

*** exports share
replace exports=exports*1000000
replace exports=0 if exports==.

* GRAPH AS A PERCENTAGE OF TRANSACTIONS 	
gen time_reg=time

collapse   tot* exports , by(time)
keep if time>=-4 & time<=4

gen x1=tot_dom + tot_gov + tot_MNC + tot_fore_nonmnc + exports
gen x2=x1-tot_gov
gen x3=x2-tot_MNC
gen x4=x3-tot_fore_nonmnc
gen x5=x4-exports

**mill dollars
forvalue i=1/5{
replace x`i'=x`i'/500/1000000
}
}
twoway rarea x1  x2 time || ///
rarea x2 x3 time || ///
rarea x3 x4 time || ///
rarea x4 x5 time || ///
area x5 time , ///
legend(order(1 "Government" 2 "MNC" 3 "Foreign (not MNC)" 4 "Exports"  5 "Domestic")) ///
 ytitle("Decomposition of Sales" , size(medlarge)) graphregion(fcolor(white)) xtitle("Event time", size(medsmall)) ylabel(  , labs(medsmall) format(%9.2fc) gstyle(shortdash)) ///
 legend(region(lwidth(none))   size(medsmall)) xtitle("Years Since First MNC Interaction", size(medlarge)) xlabel(,labsize(large))
graph export "${path0}/results/3-appendix_c/Supplementary_Figure_C2.eps", as(eps) replace

end

********************************************************************************
* Tables C2 and C3
* Heterogeneity by Number of MNC Events Experienced
********************************************************************************
prog multiple_events
display "multiple events"
quiet{
*using previous data
use $analysis_data, clear
ds D_m4-D_4
global dummies `r(varlist)'
***PREPARING FOR TRANSACTIONS	   

* SELLERS
preserve
keep ID year event past sector4 province ingre D* sector_g l_not_MNC_event_2 all_trans
rename * s_*
rename s_ID seller
rename s_year year
compress
tempfile sellers_id
save `sellers_id', replace
restore

* MNCs 
use $MNCs, clear 
keep ID 
rename ID buyer
tempfile f_buyer
save `f_buyer'
use $revec_group , clear
keep ID ZF* sop sector_g
reshape long ZF_broad_ , i(ID) j(year)
drop if year<2008
gen ZF_tot=(ZF+ZF_b>0)
bys ID: egen Zona_F= max(ZF_t)
bys ID: gen dup=cond(_N==1,0,_n)
drop if dup>1
keep ID Zona_F* sop sector_g
rename * b_*
rename b_ID buyer 
merge 1:1 buyer using `f_buyer'
keep if _m==3
drop _m
tempfile f_buyer
save `f_buyer'

*** SELLERS AND TRANSACTIONS: TOTAL SALES 
use $transactions, clear
*some transactions appear to be larger than the reported revenues of
*the seller. Those are labelled as "potentially wrong." We exclude them
drop if pot==1
drop pot
bys seller year: egen tot_trans=total(trans_)
* TRANSACTIONS TO REAL VALUE
local trans_vari "trans tot_trans"
foreach v of local trans_vari{
forvalue i=2008/2017{
replace `v'=(100*(`v')/${ipc`i'})/1000000 if year==`i'
}
}

* sellers transactions
merge m:1 seller year using `sellers_id'
drop if _m==1
drop _m
bys seller year: egen double all_trans=total(trans)
replace all_trans=s_ingre if all_trans>s_ingre & all_trans!=.
}
quiet{
*MNC buyers
merge m:1 buyer using `f_buyer'
gen fore_buyer=(_m==3)
drop if _m==2
drop _m

*saving temp file 
save temp/temp_trans.dta , replace
 
***
*** first, second, and third events
***
*first event
gen temp_year_event=year*fore_buyer
replace temp_year_event=. if temp_year_event==0
bys seller: egen event1=min(temp_year_event)
replace event1=. if event1==2008 | event1==2009
assert event1==s_event
drop temp*
gen temp_buyers=fore_buyer*(event1==year)
bys seller buyer: egen temp_trig=max(temp_buyers)
*second event
gen temp_year_event=year*fore_buyer
replace temp_year_event=. if temp_year_event==0 | temp_trig==1
bys seller: egen event2=min(temp_year_event)
replace event2=. if event2==2008 | event2==2009
assert event2>event1 | event1==. | event2==.
drop temp*
gen temp_buyers=fore_buyer*(event2==year | event1==year)
bys seller buyer: egen temp_trig=max(temp_buyers)
*third event 
gen temp_year_event=year*fore_buyer
replace temp_year_event=. if temp_year_event==0 | temp_trig==1
bys seller: egen event3=min(temp_year_event)
replace event3=. if event3==2008 | event3==2009
assert event3>event2 | event3==. | event2==.
drop temp*

collapse all_trans event*, ///
 by(seller year )
rename seller ID 
merge 1:1 ID year using $analysis_data
assert _m==3
drop _merge

***
*** regressions for three diff samples: firms receiving one event, two, or
*** three events or more 
***
global var1 "y"
global name1 "sales"
global var2 "y l k m"
global name2 "CD"

estimates clear 
local j=1

forvalue v=1/2{

*baseline 
global REG "reghdfe ${var`v'}  D_* if past!=1 , absorb(year#sector4#province IDs) vce(cluster province#sector_g)"
run_reg
addfe Never_Matched
estimates store m`j'
local j=`j'+1
local label_ "`label_' "${name`v'} base""
*only one event*
global REG "reghdfe ${var`v'}  D_* if past!=1 & event2==. & event3==., absorb(year#sector4#province IDs) vce(cluster province#sector_g)"
run_reg
addfe Never_Matched
estimates store m`j'
local j=`j'+1
local label_ "`label_' "${name`v'} 1 event""
*two events event*
global REG "reghdfe ${var`v'}  D_* if past!=1 & ((event2!=. & event3==.) | event==.)  , absorb(year#sector4#province IDs) vce(cluster province#sector_g)"
run_reg
addfe Never_Matched
estimates store m`j'
local j=`j'+1
local label_ "`label_' "${name`v'} 2 event""
*three events or more*
global REG "reghdfe ${var`v'}  D_* if past!=1 & (event3!=. | event==.), absorb(year#sector4#province IDs) vce(cluster province#sector_g)"
run_reg
addfe Never_Matched
estimates store m`j'
local j=`j'+1
local label_ "`label_' "${name`v'} >=3 event""
}
* Prepare estimates for -estout-
	estfe m*, labels(IDs "Firm FE" year#sector4#province "Year-4DSect-Prov FE")
	return list	
}	
esttab m* using ${path0}/results/3-appendix_c/Supplementary_Table_C2.tex , ///
 indicate(`r(indicate_fe)' Never_Matched)  star(* 0.10 ** 0.05 *** 0.01) varwidth(23) ///
 keep( $dummies ) cells("b(fmt(3) star )" "se(par fmt(3))" )  ///
 stats(r2_a mean_dep sd_dep n_obs n_fe n_firms , fmt(a2) labels("Adjusted R$^2$" "Mean Dep. Var. (level)" "SD Dep. Var. (level)" "\# Observations" "\# Fixed Effects" "\# Firms"))   ///
 mlabel("At least 1 event year" "1 event year only" "2 event years only" "3 event years only" "At least 1 event year" "1 event year only" "2 event years only" "3 event years only") alignment(c) ///
 replace noconstant label nodepvar collabels(none) compress type 

**** 
**** OUTCOMES BEF VS AFT FOR THE FIRST, SECOND, THIRD EVENT
**** 
quiet{
gen after1=(year>=event1)
gen after2=(year>=event2)
gen after3=(year>=event3)
gen after=.
estimates clear 
local j=1
local label_ ""
forvalue v=1/2{

*first event  
replace after=after1
global REG "reghdfe ${var`v'}  after if past!=1 , absorb(year#sector4#province IDs) vce(cluster province#sector_g)"
run_reg
addfe Never_Matched
estimates store m`j'
local j=`j'+1
local label_ "`label_' "${name`v'} base""
*second event*
replace after=after2
global REG "reghdfe ${var`v'}  after if past!=1 , absorb(year#sector4#province IDs) vce(cluster province#sector_g)"
run_reg
addfe Never_Matched
estimates store m`j'
local j=`j'+1
local label_ "`label_' "${name`v'} 2nd event""
*third event*
replace after=after3
global REG "reghdfe ${var`v'}  after if past!=1  , absorb(year#sector4#province IDs) vce(cluster province#sector_g)"
run_reg
addfe Never_Matched
estimates store m`j'
local j=`j'+1
local label_ "`label_' "${name`v'} 3rd event""
}
* Prepare estimates for -estout-
	estfe m*, labels(IDs "Firm FE" year#sector4#province "Year-4DSect-Prov FE")
	return list	
}	
esttab m* using ${path0}/results/3-appendix_c/Supplementary_Table_C3.tex , ///
 indicate(`r(indicate_fe)' Never_Matched)  star(* 0.10 ** 0.05 *** 0.01) varwidth(23) ///
 keep( after ) cells("b(fmt(3) star )" "se(par fmt(3))" )  ///
 stats(r2_a mean_dep sd_dep n_obs n_fe n_firms , fmt(a2) labels("Adjusted R$^2$" "Mean Dep. Var. (level)" "SD Dep. Var. (level)" "\# Observations" "\# Fixed Effects" "\# Firms"))   ///
 mlabel("After first event year" "After second event year" "After third event year" "After first event year" "After second event year" "After third event year") alignment(c) replace noconstant label nodepvar collabels(none) compress type 
end


********************************************************************************
* Table C4
* Heterogeneity in TFP Gains Based on the Characteristics of the First MNC Buyer
********************************************************************************
prog tfp_mnc_char_hete
quiet{
estimates clear
*size as total sales
use $revec_group , clear 
keep if foreign_list==1 
keep ID ingresosir* imports* exports* trabaj_*
local variables_ "ingresosir imports exports"
*real 
forvalue year =2008/2017{
foreach vari of local variables_{
replace `vari'_`year'= (100*(`vari'_`year')/${ipc`year'})/1000
format `vari'_`year' %16.0g
}
}
reshape long ingresosir_ imports_ exports_ trabaj_ , i(ID) j(year)
drop if trabaj_==.
bys ID: egen b_entry=min(year)
drop if year<2008 
replace imports=0 if imports==.
replace exports=0 if exports==.
*approximation to US dollars 
replace ingresos=ingresos/500
replace ingresos=exports if exports> ingresos 
collapse b_ave_sales=ingresos b_ave_imports=imports b_ave_trab=trabaj_ b_entry , by(ID)

rename ID buyer 
tempfile size_sales
save `size_sales', replace

*size as transactions purchases 
use $transactions , clear
drop if pot==1 
*real 
forvalue year =2008/2017{
replace trans= (100*trans/${ipc`year'})/1000/500 if year==`year'
}
collapse (sum) trans , by(buyer year)
collapse b_ave_purchases=trans , by(buyer)
tempfile size_trans
save `size_trans', replace

*merging with analysis data
use $analysis_data , clear

*merging sales of buyer 
merge m:1 buyer using `size_sales'
drop if _m==2
drop _m
*merging transactions of buyer 
merge m:1 buyer using `size_trans'
drop if _m==2
drop _m

tab event 
summ b_ave* if event!=.
gen share_purchases=b_ave_purchases / (b_ave_purchases+b_ave_imports) 

foreach v in b_ave_trab b_ave_purchases share_purchases b_ave_sales  {
*above/below median
global med_unit "buyer"
global med_var "`v'"
quiet ab_med
}
***
*** regressions 
*** 

foreach v in b_ave_trab b_ave_sales b_ave_purchases share_purchases  {
forvalue i=0/1{
global REG "reghdfe y k l m  D_* if event!=. & ab_med_`v'==`i' , absorb(year#sector2#province  IDs) vce(robust)"
run_reg
estimates store `v'_`i'
}
}
}
esttab b_ave_trab* b_ave_sales* b_ave_purchases* share_purchases*  ///
 using ${path0}/results/3-appendix_c/Supplementary_Table_C4.tex , ///
  star(* 0.10 ** 0.05 *** 0.01) keep(${dummies}) cells("b(fmt(2) star )" "se(par fmt(2) )" ) ///
 stats(r2_a mean_dep sd_dep n_obs n_fe n_firms , fmt(a2) labels("Adjusted R$^2$" "Mean Dep. Var. (level)" "SD Dep. Var. (level)" "\# Observations" "\# Fixed Effects" "\# Firms" ))   ///
 mlabel("Below Med. Workers" "Above Med. Workers" "Below Med. Sales" "Above Med. Sales" "Below Med. Purch." "Above Med. Purch." "Below Med. Sh. Purch." "Above Med. Sh. Purch.") alignment(c) ///
 replace noconstant label nodepvar collabels(none) compress type
end

********************************************************************************
* Table C5
* Heterogeneity in TFP Gains Based on the HQ Country of the First MNC Buyer
********************************************************************************
prog het_mnc_countries

quiet{
*loading data	
use $analysis_data , clear

*dummies to show in tables  	
ds D_m4-D_4
global dummies `r(varlist)'

*merging countries of MNCs
merge m:1 buyer using $mnc_countries
summ _m if event!=.
*assert `r(mean)'==3
drop _m

*saving data
save temp/country_info.dta , replace
}
**********************************
*   REG BY CONTINENT OF MNC 
**********************************
quiet{
*use country info 
use temp/country_info.dta , clear

*continents
global usa_can "US CA"
global europe "ES CH GB DE SV FR IT LU NL DK IE BE SE NO AT HU RS GR BG"
global asia_aus "CN JP ID KR SG HK IN AU"
global latin "PA MX CO NI CL VE GT BR PE CW HN TT BZ EC BM KY"
gen continent=""

foreach v of global usa_can {
replace continent="USA-Canada" if country=="`v'"
} 
foreach v of global europe {
replace continent="Europe" if country=="`v'"
}
foreach v of global asia_aus {
replace continent="Asia-Aust" if country=="`v'"
}
foreach v of global latin {
replace continent="Latinoamerica" if country=="`v'"
}

encode continent, gen(cont)
global categ "cont"

levelsof ${categ}, local(levels)
display "`levels'" 
local lbe: value label ${categ}
display "`lbe'"

foreach l of local levels{
local f`l': label `lbe' `l'
}
local models ""
forvalue l=1/4{
global REG "reghdfe y k l m  D_* if ${categ}==`l' , absorb(year#sector2#province  IDs) vce(robust)"
run_reg
estimates store ${categ}_`l' 
}
}  
**********************************
*   MANAGEMENT PRACTICES 
**********************************
quiet{
*management practices data 
import excel  $management , sheet("Sheet1") firstrow clear

*compute weighted average per country 
gen w=.
replace w=50 if employment=="A) 50 to 100"
replace w=100 if employment=="B) 101 to 250"
replace w=250 if employment=="C) 251 to 500"
replace w=500 if employment=="D) 501 to 1000"
replace w=1000 if employment=="E) 1000+"

*increasing weight when observed several times
replace w=w*N

*weighted average by country
collapse (mean) management [aweight=w], by(country)

*saving 
tempfile manag
save `manag',replace

*loading country codes
import excel $country_codes , sheet("Hoja1") firstrow clear
*sudan appears repeated. It does not matter because we do not use it
drop if country=="Sudan"
*renaming before the merge
replace country="United States" if country=="United States of America"
replace country="Vietnam" if country=="Viet Nam"
replace country="Tanzania" if country=="Tanzania, United Republic of"
replace country="Republic of Ireland" if country=="Ireland"
replace country="Great Britain" if country=="United Kingdom"
rename country_name country

*merging with management data 
merge 1:1 country using `manag'
drop if country=="Northern Ireland"
drop if _m==1
summ _m 
assert `r(mean)'==3
rename country country_name
rename code2 country
keep country* manage
*saving 
tempfile manag
save `manag',replace

*use events data
use temp/country_info.dta , clear
merge m:1 country using `manag'
tab country if _m==1
tab country if _m==3
keep if _m==3
drop if year==.

*regressions variable: country with above median management score
*above/below median 
global med_unit "country"
global med_var "management"
quiet ab_med

*regressions 
forvalue v=0/1{
global REG "reghdfe y k l m  D_* if event!=. & ab_med==`v' , absorb(year#sector2#province  IDs) vce(robust)"
run_reg
estimates store manag_above_`v'
}
}
collapse management ab_med, by(country*)
sort management
list
**********************************
* Reg above/below gdp pc ppp
**********************************
quiet{
*loading country codes
import excel $country_codes , sheet("Hoja1") firstrow clear
 bys code2: gen dup=cond(_N==1,0,_n)
 drop if dup>1
rename code2 country 
keep country*
tempfile country_name
save `country_name', replace

*use country info 
use temp/country_info.dta , clear
merge m:1 country using `country_name'
keep if _m==3 & country!="VI"

*above/below median
global med_unit "country"
global med_var "gdp"
quiet ab_med


*regressions 
forvalue v=0/1{
global REG "reghdfe y k l m  D_* if event!=. & ab_med==`v' , absorb(year#sector2#province  IDs) vce(robust)"
run_reg
estimates store gdp_above_`v'
}
}

**********************************
* REGRESSION TABLE
**********************************
esttab ${categ}_4 ${categ}_2 ${categ}_3 ${categ}_1 ///
gdp_above_0 gdp_above_1 manag_above_0 manag_above_1 ///
 using ${path0}/results/3-appendix_c/Supplementary_Table_C5.tex , ///
  star(* 0.10 ** 0.05 *** 0.01) keep(${dummies}) cells("b(fmt(2) star )" "se(par fmt(2) )" ) ///
 stats(r2_a mean_dep sd_dep n_obs n_fe n_firms , fmt(a2) labels("Adjusted R$^2$" "Mean Dep. Var. (level)" "SD Dep. Var. (level)" "\# Observations" "\# Fixed Effects" "\# Firms" ))   ///
 mlabel("USA-Canada" "Europe" "Asia-Pacific" "LAC" "Below Median GDPpcPPP}" "Above Median GDPpcPPP" "Below Median Manag" "Above Median Manag") ///
 alignment(c)  replace noconstant label nodepvar collabels(none) compress type
collapse gdp ab_med, by(country* )
sort gdp
list
end

********************************************************************************
*run regressiona and stores relevant info
********************************************************************************
prog run_reg
capture drop resid
capture drop dep

*regression 
$REG resid(resid) 

* add number of observations 
scalar n_obs=e(N)
estadd scalar n_obs

* add distinct number of firms
distinct IDs if resid!=.
scalar n_firms=r(ndistinct)
estadd scalar n_firms 

* add mean and sd of dependent variable 
global dep_var "`e(depvar)'"
if "$quant" == ""{
*in US dollars
gen dep=exp( $dep_var )/500
}
if "$quant" == "1" {
gen dep=exp( $dep_var )
}
global quant ""
quietly summ dep if resid!=.
scalar mean_dep=r(mean)
scalar sd_dep=r(sd)
estadd scalar mean_dep
estadd scalar sd_dep

*add number of fixed effects 
scalar n_fe=e(df_a_initial)
estadd scalar n_fe

capture drop resid
capture drop dep
end

********************************************************************************
*above/below median
********************************************************************************
prog ab_med 

egen tag=tag($med_unit )
gen temp_tag=tag*$med_var 
replace temp_tag=. if temp_tag==0
summ temp_tag, det
local med `r(p50)'
gen temp_ab_med = ($med_var >=`r(p50)') if tag!=0
replace temp_ab_med=. if $med_var ==.
bys $med_unit : egen ab_med_$med_var = max(temp_ab_med)
capture drop temp* tag  

end

quiet{
***  PRODUCER PRICE INDEX  (yearly average) 2013 ==100   
global ipc2005 = 50.99*100/92.88
global ipc2006 = 57.96*100/92.88
global ipc2007 = 63.38*100/92.88
global ipc2008 = 75.43*100/92.88
global ipc2009 = 77.53*100/92.88
global ipc2010 = 81.15*100/92.88
global ipc2011 = 85.57*100/92.88
global ipc2012 = 89.78*100/92.88
global ipc2013 = 92.88*100/92.88
global ipc2014 = 98.99*100/92.88
global ipc2015 = 98.93*100/92.88
global ipc2016 = 100.08*100/92.88
global ipc2017 = 100.77*100/92.88
}

********************************************************************************
*PROGRAM TO ADD FIXED EFFECTS INDICATOR TO LATEX TABLES           
********************************************************************************
program addfe
    if `"`0'"'=="" local 0 _fe
    tempname b V
    mat `b' = 0
    mat coln `b' = `0'
    mat `b' = e(b), `b'
    mat `V' = e(V)
    mat `V' = (`V', J(rowsof(`V'), 1, 0)) \ (J(1, colsof(`V'), 0), 0)
    erepost b=`b' V=`V', rename
end

******************************************************************************** 
***  EXECUTE THE PROGRAM main   
********************************************************************************
main

******************************************************************************** 
***  TIMER
********************************************************************************
timer off 1
*time in seconds
timer list 1
*time in minutes
local time_ =round(`r(t1)'/60)
display "the code takes `time_' minutes in total"

log close
translate `logfile'.smcl `logfile'.txt , replace linesize(250)
capture erase `logfile'.smcl
